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Abstract 

The  paper  concerns  the  modeling  of  turbulent  liquid  oxygen/hydrogen  spray  combus- 
tion for  elevated  subcritical  pressure  and  cryogenic  inlet  temperature  conditions.  Various 
approaches  are  outlined  and  discussed  that  concern  current  and  future  models  for  turbu- 
lent two-phase  flows  as  well  as  models  to  include  detailed  chemical  reactions.  The  pres- 
ence of  the  liquid  phase  complicates  the  situation  since  the  turbulence  and  the  chemical 
reactions  not  only  interact  with  each  other  but  also  with  the  spray  processes.  After  the 
presentation  and  discussion  of  general  approaches,  the  combustion  in  a single  injector 
combustion  chamber  is  modeled  where  experimental  data  are  available  for  gas  phase  tem- 
perature and  both  droplet  size  and  velocities.  The  model  uses  an  Eulerian-Lagrangian 
formulation  for  the  gas  and  the  liquid  phase,  respectively.  Detailed  models  for  droplet 
heating  and  vaporization  in  a convective  flow  field  are  employed,  and  detailed  gas  phase 
reactions  are  accounted  for  through  use  of  a flamelet  model  for  turbulent  spray  combus- 
tion. The  results  show  a very  good  agreement  between  experimental  and  computational 
spray  characteristics.  The  computed  gas  phase  temperature  lies  above  the  experimental 
values  which  is  associated  with  CARS  single  shot  measurements  and  incomplete  data  for 
the  initial  conditions  of  the  combustion  process. 

1 Introduction 

An  improved  understanding  of  the  physical  and  chemical  processes  occurring  in  liquid  rocket 
engines  is  required  to  ensure  the  stability,  reliability,  and  efficiency  of  their  performance.  The 
gaseous  hydrogen  and  the  liquid  oxygen  (LOX)  are  injected  at  cryogenic  inlet  temperatures, 
and  the  turbulent  combustion  occurs  in  both  the  sub-  and  supercritical  domain.  Therefore, 
the  models  for  the  processes  in  systems  such  as  the  Ariane  V or  the  Space  Shuttle  main 
engine  are  very  complex.  The  paper  concerns  principal  approaches  to  the  modeling  of  these 
processes  and  discusses  the  state  of  the  art  as  well  as  potential  approaches  in  future. 

Finally,  the  combustion  process  in  a single  injection  combustion  chamber  is  investigated 
for  elevated  pressure  and  cryogenic  inlet  temperatures.  Experimental  data  are  available  for 
an  elevated  pressure  situation  of  5 bar,  and  the  experimental  data  are  published  by  Sender  et 
al  [1],  The  combustion  in  liquid  rocket  propulsion  typically  occurs  in  the  flamelet  regime  of 
turbulent  combustion  [2]  which  enables  the  use  of  the  flamelet  model  for  turbulent  spray  dif- 
fusion flames  [3,  4].  Structures  of  laminar  hydrogen/oxygen  flames  are  precalculated  at  5 bar 
and  for  cryogenic  inlet  temperature  of  100  K for  hydrogen  [5,  6].  They  are  incorporated  into 
the  turbulent  spray  combustion  model  which  allows  the  computation  of  all  chemical  species 
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that  are  included  in  the  detailed  chemical  reaction  mechanism  for  the  hydrogen/oxvgen 
system  that  includes  8 reactive  species  and  38  chemical  reactions  [7]. 


2 Turbulence  Modeling  of  Two-Phase  Flows 

Two-phase  flows  are  characterized  by  the  separation  of  the  turbulent  flow  field  and  the  strong 
interaction  between  the  spray  and  the  gas  through  processes  such  as  turbulent  mixing,  droplet 
heating,  vaporization,  and  motion  as  well  as  coupling  with  possible  occurance  of  chemical 
reactions.  In  practical  systems,  most  often  the  turbulent  gas  flow  is  modeled  using  a fc  — e 
turbulence  model  where  additional  terms  are  added  to  account  for  the  interaction  with  the 
liquid.  This  formulation  is  done  within  an  Eulerian  Lagrangian  formulation  of  governing 
equations  for  a dilute  spray  within  either  a discrete  droplet  or  a continuous  droplet  model, 
c.f.  Faeth  [8].  An  example  of  the  mathematical  equations  within  a discrete  droplet  model  is 
given  through  [3] 
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where  the  source  terms  are  given  in  Tab.  1.  Here  a two  dimensional  axisymmetric  configura- 
tion is  considered  where  x and  r denote  the  axial  and  radial  directions,  respectively.  Favre 
averaged  values  are  used. 
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Table  1:  Source  terms  for  Eq.  1 [3]. 
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Table  2:  Source  terms  for  Eq.  2 [12]. 


Even  though  this  model  is  used  in  most  technical  simulations,  it  has  major  drawbacks  such 
as  the  well  known  problem  of  the  k — e model  in  strongly  recirculating  flows  and  its  failure  of 
modeling  anisotropic  turbulence  which  is  common  to  all  models  using  the  turbulence  energy 
k.  Moreover,  the  model  is  based  on  a gradient  diffusion  assumption  which  is  not  appropriate 
to  predict  the  effect  of  counter  gradient  diffusion  which  has  been  observed  in  turbulent  gas 
flows. 

More  advanced  models  that  are  under  investigation  in  gas  phase  flows  are  currently 
extended  to  two-phase  flows.  These  are  the  LES  - Large  Eddy  Simulation,  the  Reynolds 
stress  models,  and  DNS  - Direct  Numerical  Simulation.  There  is  a special  presentation  [9]  at 
the  present  meeting  that  deals  with  LES,  and  therefore  this  approach  is  not  discussed  within 
the  frame  of  this  presentation.  The  method  of  DNS  is  currently  applied  to  small  problems 
- small  in  the  sense  of  physical  dimensions  of  the  system  investigated  [10]  since  the  entire 
range  of  length  scales  needs  to  be  resolved  which  restricts  the  dimension  of  the  problem. 
However,  the  method  is  very  interesting  when  well  defined  submodels  are  to  be  investigated. 

The  Reynolds  stress  models  are  typically  used  when  non-isotropic  turbulence  occurs  and 
if  counter  gradient  diffusion  is  present.  Thejnodel  does  not  use  the  k equation  but  transport 
equations  for  the  Reynolds  stress  terms  u'-Uj.  Additional  terms  appear  that  describe  the 
coupling  between  the  velocity  fluctuations  and  the  spray  source  terms  that  stem  from  droplet 
heating,  vaporization,  and  motion.  The  Reynolds  stress  equations  yield 
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and  the  source  terms  are  given  in  Tab.  2: 
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For  completeness  and  further  use  in  the  next  section,  the  element  mass  fraction.  Zj , is 
also  shown  in  the  Tab.  2.  The  terms  I to  IV  in  Eq.  2 may  be  closed  using  the  standard 
models  for  pure  gas  phase  equations  whereas  ther  terms  V describe  the  interaction  between 
the  spray  and  the  variable  under  consideration.  These  terms  may  be  formulated  according 
to  the  derivations  in  Ref.  [3,  4],  and  they  are  currently  under  investigation. 

Another  approach  is  that  of  Burlaka  mid  Borghi  [11]  that  describes  the  spray  surface. 
with  a transport  equation 
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The  destruction  term  includes  the  vaporization  process.  This  model  is  motivated  by  the 
surface  transport  equation  in  turbulent  premixed  combustion,  and  the  description  of  the 
vaporization  process  is  still  under  development. 


3 Models  for  Chemical  Reactions  in  Turbulent  Reactive  Two- 
Phase  Flows 


The  chemical  reactions  in  a convective  flow  field  may  be  described  mathematically  through 
the  conservation  equation  of  the  mass  fractions,  V;,  of  chemical  species  i : 
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where  ibi  is  the  molar  chemical  reaction  rate  and  Lv  is  the  source  therm  due  to  liquid 
vaporization.  Sn  is  the  Dirac  Delta  function  where  L denotes  the  species  in  liquid  phase. 
The  system  of  chemical  reactions  is  written  as 
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and  the  consumption/production  rate  of  a chemical  species  yields 
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where  ivj  is  the  reaction  rate  of  a reaction  step,  j,  using  the  modified  Arrhenius  expression 
for  the  reaction  constant 
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for  N chemical  reaction  steps. 

If  the  chemical  reactions  are  to  be  modeled  using  a mixture  fraction,  the  source  term  due 
to  chemical  reactions  disappears.  The  mixture  fraction  may  be  based  on  the  elements,  d 
present  in  the  system: 
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where  i = N,  O,  H,  and  fitJ  is  the  mass  of  element  i in  molecule  j.  The  mixture  fraction  may 
also  be  based  on  enthalpy  [13] 
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If  the  Lewis  numbers  of  all  species  were  unity,  then  there  is  no  difference  in  the  above  defini- 
tions. However,  for  the  hydrogen/oxygen  system,  this  is  not  true,  and  some  of  the  definitions 
show  unreasonable  values  [13].  In  the  remainder  of  the  paper,  the  mixture  fractions  is  based 
on  the  element  H which  is  the  best  choice  in  terms  of  monotonicity  with  physical  space  for 
the  hydrogen/oxygen  system. 

Taking  these  definitions  for  the  instantaneous  values  of  the  chemical  reaction  rate,  there 
axe  various  approaches  to  close  the  problem  in  turbulent  reactive  flows.  The  directest  ap- 
proach is  the  direct  closure  of  the  source  term  of  chemical  reactions  through  use  of  a joint 
probability  density  function  (PDF),  P,  leading  to  the  following  averaged  chemical  reaction 
rate: 

/■]  rl  roc  roc 
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Since  the  joint  PDF  is  not  known,  most  often  it  is  factorized  and  beta  functions  are  used  to 
describe  to  single  PDFs 
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The  a and  j3  are  linear  functions  of  the  variance  of  the  variable  $ under  consideration,  and 
they  are  obtained  through  solution  of  their  transport  equations,  c.f.  Eq.  2. 

Most  often  this  approximation  is  not  too  bad  in  gas  phase  combustion.  However,  it  ap- 
prears  [10]  that  the  beta-function  is  a poor  approximation  in  regimes  where  the  vaporization 
occurs.  Therefore,  new  approaches  need  to  be  developed. 

A promising  approach  seems  to  be  the  formulation  through  a PDF  transport  equation. 
This  model  has  been  developed  by  Pope  [14]  and  Dopazo  [15]  for  gas  phase  combustion.  The 
PDF  is  modeled  through  a transport  equation,  and  the  chemical  source  terms  may  be  used 
without  the  averaging  procedure  described  above.  The  formulation  of  the  PDF  transport 
equation  for  the  PDF  of  the  mixture  fraction,  yields  for  turbulent  two-phase  flows 


dt 
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where  the  third  term  on  the  LHS  is  usually  modeled  using 

{P<u=m=^m  (is> 

This  approach  has  not  yet  been  investigated  for  turbulent  two-phase  flows  to  the  author’s 
knowledge,  but  it  appears  promising.  The  diffusion  term  needs  to  modeled  using  particle 
mixing  models,  and  it  is  unclear  if  they  can  be  used  for  two-phase  flows  without  modification. 
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The  major  advantage  of  this  approach  is  that  the  PDF  tranport  equation  needs  no  infor- 
mation on  fluctuations  of  the  droplet  vaporization,  and  it  enables  the  evaluation  of  a PDF 
that  may  possibly  be  approximated  through  an  extension  of  well  known  standard  PDFs. 
For  reactive  flows,  the  approach  is  beneficial  because  it  does  not  require  averaged  chemical 
reaction  terms  the  instantaneous  values  are  sufficent. 

The  chemistry  itself  may  be  included  through  use  of  detailed  reaction  schemes  [7].  reduced 
reaction  mechanisms  [16],  or  tabulated  systems  such  as  the  ILDM  approach  [17].  There  is 
a number  of  reduced  chemcial  reaction  schemes  for  the  hvdrogen/oxygen  system,  and  they 
range  from  4-step  to  2-step  schemes  [16].  The  chemcial  reactions  in  high  pressure  systems 
typically  occur  with  small  time  scales  so  that  the  use  of  a lower  number  of  reaction  steps  is 
more  justified  as  pressure  increases,  and  for  very  high  pressure  even  a one-step  mechanism 
performs  well. 

Another  type  of  model  to  include  detailed  chemical  reactions  is  the  flamelet  model  for 
turbulent  combustion  [18,  19,  20]  which  has  been  extended  for  turbulent  spray  diffusion 
flames  [3.  4].  The  model  considers  the  tubulent  flame  to  consist  of  an  ensemble  of  laminar 
flamelets  that  are  stretched  through  the  turbulence  of  the  flow  field,  c.f.  Fig.  1. 
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Figure  1:  Schematic  of  the  flamelet  model  for  turbulent  spray  combustion. 

The  laminar  flamelet  depends  on  the  mixture  fraction.  the  scalar  dissipation  rate,  y, 
the  initial  droplet  size  and  velocity,  Rio  and  vio , resp..  and  the  equivalence  ratio,  Er . The 
Favre  averaged  values  of  a variable,  <f>.  then  yields  [4] 
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The  joint  probability  density  function  is  factorized,  and  the  one-dimensional  PDFs  need  to 
be  determined  [4].  Currently,  there  are  libraries  available  for  the  liquid  oxygen/ hydrogen 
system  [6]  and  for  the  methanol/air  [21]  and  ethanol/air  [22]  systems. 

In  the  following  section,  the  flamelet  model  for  turbulent  spray  diffusion  flames  is  used 
together  with  the  model  described  through  the  equations  given  in  Tab.  1 to  simulate  the  com- 
bustion process  in  the  micro  combustion  chamber  M3  of  the  DLR  in  Lampoldshausen  [23]. 
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4 Simulation  of  the  Combustion  Process  in  the  Micro  Com- 
bustion Chamber  M3 

4.1  Model 

A dilute  LOX  spray  is  considered  that  is  injected  into  a turbulent  gaseous  hydrogen  stream 
where  the  inlet  temperatures  are  cryogenic.  The  system  pressure  is  5 bar.  The  model  in- 
cludes an  Eulerian  description  of  the  gas  phase  and  Lagrangian  equations  for  the  dilute 
spray.  The  k — e turbulence  model  is  employed  where  additional  terms  account  for  the  spray 
interaction  [3],  c.f.  Tab.  1.  The  chemical  reactions  are  described  through  a flamelet  model 
for  turbulent  spray  diffusion  flames  [3,  4]  outlined  in  the  previous  section.  The  conserva- 
tion equations  for  the  mixture  fraction  and  its  variance  also  account  for  mass  gain  through 
vaporization  of  the  liquid  [3]. 

Convective  heating  and  vaporization  is  described  through  a model  developed  by  Abram- 
zon  and  Sirignano  [24].  The  equation  for  droplet  motion  accounts  for  turbulence  effects 
through  a Gaussian  distribution  for  turbulent  fluctuations  [25].  The  spray  distribution  is 
described  through  the  discrete  droplet  model  [23]. 

For  the  computation  of  the  gas  phase  characteristics  in  the  cryogenic  temperature  regime, 
data  from  JSME  tables  [26]  for  pressures  up  to  200  bar  and  temperatures  between  80  and 
300  K are  used.  Moreover,  the  pressure  (and  temperature)  dependence  of  the  vaporization 
rate  and  of  the  bina^  equilibrium  composition  at  the  liquid/gas  interface  [27]  is  included. 
For  the  pressure  considered  here,  the  real  gas  effects  in  the  gas  phase  are  negligable  [28] . 

4.2  Results  and  Discussion 


Figure  2:  Measured  droplet  velocities  in  the  range  between  72  and  140  mm  from  their 
injector  [1]. 

The  experimental  data  [1,  23]  include  gas  temperature  as  well  as  droplet  velocites  and 
size.  Figure  2 displays  the  droplet  velocities  obtained  from  experiment  using  PIV.  The  area 
displayed  ranges  from  72  mm  to  140  mm  of  the  combustion  chamber  which  is  the  second  half 
of  the  entire  single-injector  micro  combustion  chamber  M3  of  the  DLR  in  Lampoldshausen, 
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Figure  3:  Initial  droplet  rate  and  distribution  for  9 different  size  classes  evaluated  from 
experimental  data  [1]  at  x = 72  mm  from  the  injector  [23]. 
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Figure  4:  Contour  plot  of  the  gas  temperature  [23]. 


Germany.  In  the  first  part  of  the  combustor,  the  spray  is  rather  thick,  and  the  model  is  not 
suitable.  The  experimental  data  are  used  to  produce  initial  profiles  for  the  simulation  as 
well  as  for  comparison  further  downstream  in  the  combustion  chamber. 

Figure  3 shows  the  generated  initial  droplet  rate  at  the  axial  position  72  mm.  The  depth 
of  field  in  the  experiment  is  0.5  mm  which  has  been  used  to  calculate  the  droplet  rate.  The 
experiment  typically  does  not  comprise  more  than  about  15  % of  the  entire  liquid.  Sometimes 
in  the  literature,  this  value  is  extrapolated  to  100%  and  there  are  several  methods  to  do  so. 
However,  we  preferred  to  stick  to  the  experimental  data  that  we  received  since  there  is  no 
reliable  information  about  how  much  of  the  liquid  flux  exactly  needs  to  be  compensated  for 
which  causes  a second  ambiguity  in  the  procedure. 
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Figure  5:  Radial  profiles  of  the  gas  temperature  for  three  different  axial  positions:  Symbols 
present  experimental  data  and  lines  show  computational  results  [23]. 
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Figure  6:  Radial  profiles  of  the  gas  temperature  and  species  mass  fractions  at  x = 

127  mm  [23]. 

Also,  there  are  no  gas  temperature  measurements  exactly  at  x — 72  mm,  and  the  inlet 
conditions  for  the  gas  temperature  have  been  obtained  from  interpolated  experimental  data 
at  two  surrounding  axial  positions.  Moreover,  the  experimental  data  have  been  extrapolated 
into  the  outer  regions  where  no  experimental  data  are  available. 

There  is  no  information  on  other  gas  phase  characteristics  such  as  species  concentra- 
tions from  the  experiment.  The  initial  gas  velocity  is  taken  from  the  experimental  velocity 
of  the  smallest  droplets,  and  the  species  profiles  are  estimated  from  the  experimental  gas 
temperature.  The  turbulence  quantities  k and  e are  computed  from  the  gas  phase  velocity. 

Figure  4 shows  the  computed  overall  structure  of  the  spray  flame  in  terms  of  the  computed 
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Figure  8:  Contour  plot  of  the  radial  gas  velocity  [23]. 
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gas  phase  temperature.  The  main  chemical  reaction  zone  is  located  at  the  boundary  of 
the  spray  jet  where  vaporized  oxidizer  meets  the  surrounding  hydrogen  gas  stream.  The 
maximum  flame  temperature  is  about  3000  K which  is  typical  for  these  high-pressure  flames. 

A more  detailed  investigation  of  flame  temperature  is  shown  in  Fig.  5 where  the  calculated 
and  experimental  radial  profiles  of  the  gas  temperature  are  shown  for  three  different  axial 
positions.  Symbols  show  experimental  data  and  lines  are  computational  results. 

Figure  5 shows  that  the  gas  temperature  in  the  main  reaction  zone  is  overpredicted 
by  the  computations.  This  may  be  attributed  to  both  experimental  and  computational 
uncertainties.  First  the  experiments  employed  CARS  single  shot  measurements,  and  the 
averaging  procedure  of  these  values  typically  leads  to  an  underprediction  of  the  measured 
gas  temperatures.  Moreover,  the  initial  gas  phase  profiles  needed  for  the  computations 
is  estimated  from  the  gas  temperature  profile,  and  this  procedure  is  somewhat  arbitrary. 
Another  reason  is  to  be  found  in  the  incomplete  liquid  phase  data  as  discussed  earlier.  The 
computations  use  the  uncorrected  data  from  measurements  so  that  the  vaporization  of  the 
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Figure  9:  Experimental  and  computational  results  of  the  Sauter  mean  radius  at  x = 
104  mm  [23]. 


Radial  Position  y [mm] 


Figure  10:  Experimental  and  computational  results  of  the  droplet  velocities  as  well  as  cal- 
culated gas  velocities  at  x = 104  mm  [23]. 

neglected  liquid  would  reduce  the  gas  phase  temperature.  Therefore,  the  discrepancies  of 
computed  and  measured  temperature  profiles  in  the  gas  phase  is  reasonable. 

The  flamelet  model  is  suitable  to  predict  all  species  profiles  that  are  considered  in  the 
laminar  flamelet  library.  Figure  6 shows  radical  profiles  at  x = 127  mm.  The  HO2  and  H2O2 
are  present  in  the  colder  flame  region  (due  to  their  stability)  near  the  axis  of  symmetry  where 
molecular  oxygen  in  the  gas  phase  is  present.  The  location  of  the  H2O  peak  is  shifted  to  the 
fuel  side  which  is  typical  for  gas  phase  combustion.  Figure  7 shows  the  contour  plot  of  the 
HO2  mass  fraction,  and  the  differences  compared  to  the  gas  temperature  contour  plot  are 
obvious.  The  HO2  prevails  into  the  colder  regions  of  the  combustion  chamber  favorably  in 
regions  where  gaseous  O2  is  present. 
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The  gas  phase  velocities  show  that  there  is  a region  where  recirculation  occurs.  This 
may  he  seen  in  Fig.  8 where  the  radial  gas  phase  velocity  is  displayed.  This  finding  is  in 
agreement  with  experiment  [1].  It  is  well  known  that  the  k — e turbulence  model  that  is 
used  here  is  not  the  best  choice  if  recirculation  is  present,  and  a future  study  will  include  a 
Reynolds  stress  model  for  turbulent  spray  flames. 

Figures  9 and  10  show  a comparison  of  liquid  phase  characteristics.  Figure  9 displays 
the  radial  profile  of  the  Sauter  mean  radius  at  x = 104  mm.  Symbols  mark  experimental 
results,  and  the  agreement  between  experiment  and  computation  is  excellent.  The  same  is 
true  for  the  profiles  of  axial  velocities  at  the  same  axial  position  plotted  in  Fig.  10.  Thus  it 
is  shown  that  the  present  model  is  suitable  to  correctly  predict  the  liquid  phase  properties. 

5 Summary  and  Future  Research 

The  combustion  process  in  the  micro  combustion  chamber  M3  has  been  studied  by  means  of 
numerical  computations.  The  Eulerian/Lagrangian  model  is  suitable  to  predict  the  processes 
in  the  dilute  spray  regime  of  the  combustion  chamber.  In  particular,  the  prediction  of  the 
spray  distribution  and  evolution  is  very  good. 

Calculated  gas  phase  temperatures  are  higher  than  the  experimental  values  which  is 
explained  by  various  reasons.  First,  the  experimental  data  are  probably  underpredicted 
which  is  typical  for  CARS  single  shot  measurements.  Moreover,  only  a small  portion  of  the 
total  liquid  flux  was  captured  by  the  experimental  techniques  so  that  the  energy  consumed 
by  droplet  vaporization  is  underpredicted  leading  to  too  high  values  of  gas  temperature  in 
the  computation.  Furthermore,  there  is  no  information  of  gas  phase  species  profiles  from 
experiment  leading  to  uncertainties  in  evaluating  the  initial  conditions  for  the  computations. 

Future  studies  should  include  the  dense  spray  regime,  and  an  extension  of  the  present 
model  into  that  region  will  be  developed. 

Moreover,  the  Reynolds  stress  model  for  turbulent  spray  flames  will  be  employed  to 
account  for  an  improved  simulation  of  the  recirculation  zone.  Here  the  coupling  terms  of 
velocities  and  vaporization  fluctuations  need  to  be  modeled. 

In  future,  the  pressure  will  also  be  increased  into  the  high-pressure  domain. 
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